LevelPool Subroutine

public subroutine LevelPool(time, dt, dtrk, Qin_t, Qin_tplusdt, pool)

Discharge routing through reservoir. Given Qin, find Qout and stage References:

Arguments

Type IntentOptional Attributes Name
type(DateTime), intent(in) :: time
integer, intent(in) :: dt

main time step (s)
integer, intent(in) :: dtrk

time step to run runge-kutta (s)
real(kind=float), intent(in) :: Qin_t

Input discharge at time t (m3/s)
real(kind=float), intent(in) :: Qin_tplusdt

Input discharge at time t + dt (m3/s)
type(Reservoir), intent(inout) :: pool

Variables

Type Visibility Attributes Name Initial
character(len=3), public :: Qdiverted
real(kind=float), public :: Qin_stepbegin
real(kind=float), public :: Qin_stepend
real(kind=float), public :: Qinpool_t

Actual input dicharge to the pool at time t (m3/s)
real(kind=float), public :: Qinpool_tplusdt

Actual input dicharge to the pool at time t+dt (m3/s)
character(len=3), public :: QoutDOY

doy to compute Qout

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
integer(kind=short), public :: doy

day of year
real(kind=float), public :: h

current stage (m)

integer, public :: j
integer, public :: nstep

number of sub-steps
logical, public :: qinIsRising

Source Code

SUBROUTINE LevelPool &
!
(time, dt, dtrk, Qin_t, Qin_tplusdt, pool)

!arguments with intent in:
TYPE(DateTime), INTENT(IN)     :: time
INTEGER, INTENT(IN)            :: dt !!main time step (s)
INTEGER, INTENT(IN)            :: dtrk !!time step to run runge-kutta (s)
REAL(KIND = float), INTENT(IN) :: Qin_t !!Input discharge at time t (m3/s)
REAL(KIND = float), INTENT(IN) :: Qin_tplusdt !!Input discharge at time t + dt (m3/s)

!arguments with intent inout:
TYPE(Reservoir), INTENT(INOUT) :: pool

!arguments with intent out:
!REAL(KIND = float), INTENT(OUT) :: Qout !!downstream discharge at time t + dt (m3/s)

!local declarations
INTEGER                         :: nstep !!number of sub-steps
REAL(KIND = float)              :: Qin_stepbegin, Qin_stepend !(m3/s)
REAL(KIND = float)              :: Qinpool_t !!Actual input dicharge to the pool at time t (m3/s)
REAL(KIND = float)              :: Qinpool_tplusdt !!Actual input dicharge to the pool at time t+dt (m3/s)
REAL(KIND = float)              :: h !!current stage (m) 
INTEGER                         :: j
LOGICAL                         :: qinIsRising
INTEGER (KIND = short)          :: doy !!day of year
CHARACTER (LEN = 3)             :: Qdiverted
CHARACTER (LEN = 3)             :: QoutDOY  !!doy to compute Qout 


!-------------------------------end of declaration-----------------------------

!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
!                                                                             !
!                                                                             !
!                                                            *   Qinpool_tplusdt  
!                                                            |                !
!                                                            |                !
!                                         Qin_stepend  O     |                !
!                                                      |     |                !
!                                                      |     |                !
!                                Qin_stepbegin   O     |     |                !
!                                             .  |     |     |                !
!                                           .    |     |     |                !
!                                         .      |     |     |                !
!                                       .        |     |     |                !
!                                     .          |     |     |                !
!                                    O           |     |     |                !
!                                    |           |     |     |                !
!                                    |           |     |     |                !
!                  Qinpool_t   * ____|...........|_____|_____|                !
!                                                                             !
!                                                <----->                      !           
!                                                  dtrk                       !
!                             t <-----------------------------> t + dt        !
!                                             dt                              !
!                                                                             !
!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!

!initialize variables to run runge-kutta
SELECT CASE ( TRIM(pool % typ) )
    
  CASE ( 'off' )
    !diverted discharge by weir is the input to the pool

    ! time t
    CALL TableGetValue ( valueIn = Qin_t, tab = pool % weir, keyIn = 'Qstream', &
                                 keyOut ='Qdiverted', match = 'linear', &
                                 valueOut = Qinpool_t, bound = 'extendconstant' )
    
    ! time t+dt
    CALL TableGetValue ( valueIn = Qin_tplusdt, tab = pool % weir, keyIn = 'Qstream', &
                                 keyOut ='Qdiverted', match = 'linear', &
                                 valueOut = Qinpool_tplusdt, bound = 'extendconstant' )
    
	  !if pool is full discharge can not enter anymore
    IF (pool % stage > pool % stageMax) THEN
      Qinpool_t = 0.
      Qinpool_tplusdt = 0.
    END IF
		
  CASE ('on')
    !discharge from stream enters pool
    Qinpool_t = Qin_t
    Qinpool_tplusdt = Qin_tplusdt
    
    IF ( Qinpool_tplusdt > Qinpool_t ) THEN
        qinIsRising = .TRUE.
    ELSE
        qinIsRising = .FALSE.
    END IF
	
	CASE DEFAULT
    CALL Catch ('error', 'LevelPool', 'wrong reservoir type in pool: ', &
	                      argument = TRIM(pool % name))
END SELECT

!initialize stage
h = pool % stage

!compute day of year (used to set ecological flow)
doy = DayOfYear (time, 'noleap') 

!set doy for computing Qout
QoutDOY = ToString ( pool % geometryDOY (doy) )

! manage high level 
IF ( pool % highLevel .AND. h > pool % fullReservoirLevel ) THEN
    
    IF ( pool % rising) THEN
        IF ( qInIsRising ) THEN
            SELECT CASE ( pool % qoutRule )
            CASE ( 1 ) !qout = qin
                 pool % Qout = Qinpool_tplusdt
            END SELECT
            RETURN
        END IF
    ELSE
        SELECT CASE ( pool % qoutRule )
        CASE ( 1 ) !qout = qin
            pool % Qout = Qinpool_tplusdt
        END SELECT
        RETURN
    END IF
    
END IF

!loop throughout steps
nstep = dt / dtrk

DO j = 1, nstep
   Qin_stepbegin = LinearInterp (0., REAL(dt), Qinpool_t, &
                                 Qinpool_tplusdt, 1.*(j-1)*dtrk)
   Qin_stepend = LinearInterp (0., REAL(dt), Qinpool_t, &
                               Qinpool_tplusdt, 1.*j*dtrk)
   
   !update level pool
   h = RungeKutta ( hini = h, Qin_begin = Qin_stepbegin, &
                    Qin_end = Qin_stepend, dt = dtrk, &
                    geometry = pool % geometry, &
                    rk = pool % rk, typOut = pool % typOut, &
                    hTarget = pool % stageTarget, &
                    eFlow = pool % eFlow (doy), &
                    freeFlow = pool % freeFlow, &
                    freeFlowElevation = pool % freeFlowElevation, &
                    QoutDOY = QoutDOY  )
END DO

!downstream discharge
SELECT CASE ( TRIM (pool % typ ) )
CASE ('on') !on-stream reservoir
    IF (pool % typOut == 2 .AND. h < pool % stageTarget) THEN
      pool % Qout = MIN( pool % eFlow (doy), Qin_tplusdt)
    ELSE IF ( pool % typOut == 2 .AND. h >= pool % stageTarget) THEN
        CALL TableGetValue ( valueIn = h, tab = pool % geometry, keyIn = 'h', &
                            keyOut ='Qout', match = 'linear', &
                            valueOut = pool % Qout, bound = 'extendconstant' )
        IF ( pool % eFlow (doy) > 0. ) THEN
                pool % Qout = MAX (pool % eFlow (doy), pool % Qout)
        END IF
    ELSE
         IF ( pool % stage <= pool % freeFlowElevation .AND. &
             Qin_tplusdt < pool % freeFlow ) THEN
            pool % Qout = Qin_tplusdt
         ELSE
            CALL TableGetValue ( valueIn = h, tab = pool % geometry, keyIn = 'h', &
                                 keyOut ='Qout', match = 'linear', &
                                valueOut = pool % Qout, bound = 'extendconstant' )
         END IF  
    END IF
  CASE ('off') !off-stream reservoir
    pool % Qout = Qin_tplusdt - Qinpool_tplusdt
END SELECT
  
!update stage
pool % stage = h
  
RETURN
END SUBROUTINE LevelPool